Analysis of tropical forest inventory data
In this section, we aimed at exploring the intraspecific variability of tree growth in three hyperdiverse mature forests, encompassing contrasting climates and disturbance regimes across the Tropics. First, for each dataset, we estimated the intraspecific variance of annualized growth using a hierarchical bayesian model. Then, we tested if spatial structure in individual mean growth could be detected. Finally, we tested if intraspecific variability was higher than interspecific variability in individual mean growth.
Note that we here used data on tree growth, which is one component of individual fitness and an imperfect proxy of tree competitive ability, and our analyses could be further tested with other traits in the future.
1 Datasets
1.1 Paracou
The Paracou forest is located in French Guiana. It is one of the best-studied lowland tropical forests in the Guiana Shield region. It belongs to the Caesalpiniaceae facies and has amongst the highest alpha-diversity in the Guiana shield with 150-200 tree species per hectare in inventories of trees with diameter at breast height (DBH) \(\ge\) 10 cm. The Guiana Shield is characterized by Pre-Cambrian granitic and metamorphic geological formations, highly eroded. It is associated with gently undulating landscapes and a very dense hydrographic system. Paracou forest lies in a hilly area, on a formation called “série Armina” characterized by schists and sandstones and locally crossed by veins of pegmatite, aplite and quartz. The topography of the site consists of small hills separated by narrow (< 5 m wide) sandy waterbeds. The altitude varies from 5 to about 45 m above sea level Hérault and Piponiot (2018). The mean annual temperature is 26 °C and winds are generally weak. There is a well marked dry season (from mid-August to mid-November) and a long rain season with a short drier period between March and April (mean annual rainfall of 3,041 mm). Different study programs have been led at the Paracou site (https://paracou.cirad.fr/). Here, we used data from a disturbance experiment, where 12 square plots of 6.25 ha were exposed to four different logging intensities between 1986 and 1988, as well as four square plots of 6.25 ha which were set up for biodiversity monitoring in 1990-1992 (plots 1-15). Since then, cartesian coordinates, DBH, species identity and survival of each tree with a DBH \(\ge\) 10 cm have been collected every one or two years, during the dry season (mid-August to mid-November). The studied 93.75-ha area harbours 91284 measured trees, among which 0 belong to an undetermined species and others belong to -1 determined species, 254 genera and 64 families.
We removed all measures of the years before the perturbations were performed and the biodiversity plots were added (1985-1991), thus obtaining a dataset containing 91029 measured trees, among which 0 belong to an undetermined species and others belong to -1 determined species, 254 genera and 64 families.
Figure 1.1 shows the location of the used trees in the Paracou site.
Figure S1.1: Plots at the Paracou site. Each point is a tree.
1.2 Uppangala
The Uppangala Permanent Sample Plot (UPSP) is located in South-East Asia, in the Western Ghats of India, and was established in 1989 by the French Institute of Pondicherry in the Kadamakal Reserve Forest, in the Pushpagiri Wildlife Sanctuary, in Karnataka state, India (Pélissier et al. 2011). It is a low altitude (500-600 m) wet evergreen monsoon Dipterocarp forest (Bec et al. 2015). This forest is considered as one of the rare undisturbed tropical forests in the world (Pascal and Pelissier 1996). The studied area of 5.07 ha is quite steep, with a mean slope angle of about 30–35°. The plots consist in five North–South oriented transects that are 20-m wide, 180- to 370-m long, and 100-m apart center to center, in addition to three rectangular plots which overlap the transects. The transects gather data from 1990 to 2011 and the rectangular plots from 1993 to 2011. The trees with GBH (Girth at Breast Height) \(\ge\) 30 cm (equivalent to ca. 9.5 cm DBH) were measured every 3 to 5 years. The original dataset contains measurements of 3870 trees of 102 species (including 2 morphospecies),0 genera and 0 families.
We removed the census dates which were not common for all plots (1990-1992).
Thus, we obtained a dataset containing 3870 trees of 102 species (including 2 morphospecies),0 genera and 0 families.
Figure 1.2 shows the location of the trees in the Uppangala site.
Figure S1.2: Plots at Uppangala site. Each point is a tree.
1.3 BCI
The Barro Colorado Island site is located in central America, in Panama, covered by a lowland tropical moist forest. The zone became an island after a valley was flooded in order to build the Panama Canal, in 1913. It is nowadays considered as the most intensively studied tropical forest in the world. The studied site is a 50 ha plot (500 \(\times\) 1,000 m). It has an elevation of 120 m and is quite flat (most slopes are gentler than 10°). Complete censuses of all trees with DBH \(\ge\) 1 cm have been performed every 5 years since 1980. The dataset contains measurements of 423617 trees of 328 tree species, 195 genera and 63 families.
When only taking into account trees with DBH \(\ge\) 10 cm for consistency with the other datasets, it contains measurements of 37224 trees of 255 tree species, 167 genera and 59 families. The dataset is available in (Condit et al. 2019).
Figure 1.3 shows the location of the selected trees in the BCI site.
Figure S1.3: The 50 ha plot of the Barro Colorado Island site. Each point is a tree.
2 Growth and mean growth estimation
For all three datasets, we computed annualized growth between two censuses as the difference of DBH between two consecutive censuses, divided by the time period between those two censuses. We then removed estimates with a growth < -2 mm/year or > 100 mm/year.
For the Paracou site, after removing 105836 growth estimates for 12520 trees of indeterminate species, we obtained a dataset with 931389 growth estimates for 65769 trees of 614 species.
For the Uppangala site, we obtained a dataset with 57921 growth estimates for 3725 trees of 102 species.
For the BCI site, we obtained a dataset with 167618 growth estimates for 30386 trees of 244 species.
For all three datasets, we finally computed the mean annual growth for each individual tree as the difference of DBH between the first and the last time a tree was measured, divided by the period between those two measurements.
Figure S2.1: Abundance diagrams for the three study sites. For the Paracou site, only individuals that were identified to the species level are accounted for. The three communities include few dominant species and many rare species (long asymptote).
3 Estimation of IV
To quantify the relative importance of intraspecific and interspecific variability for each site, we built a model including a species random effect and an individual random effect on the intercept. In order to compare to classic growth models, and estimate an intraspecific variability within diameter classes, we also added an intercept and a diameter fixed effect.
\(ln(G_{ijt} + 2) = [\beta_0 + b_{0j} + d_{0i}] + \beta_1 \times ln(D_{ijt}) + \epsilon_{ijt},\)
Priors
\(\beta_0 \sim \mathcal{N}(mean=0, var = 1), iid\)
\(\beta_1 \sim \mathcal{N}(mean=0, var = 1), iid\)
\(b_{0j} \sim \mathcal{N}(mean=0, var=V_bj), iid\)
\(d_{0j} \sim \mathcal{N}(mean=0, var=V_bi), iid\)
\(\epsilon_{ijt} \sim \mathcal{N}(mean=0, var=V), iid\)
Hyperpriors
\(V_{bj} \sim t(3, 0, 2.5), iid\)
\(V_{bi} \sim t(3, 0, 2.5), iid\)
\(V \sim t(3, 0, 2.5), iid\)
Where \(i\) is the individual and \(j\) the species. The model included two fixed effects (\(\beta_0\) and \(\beta_1\)), a random species effect \(b_{0j}\) on the intercept and a random individual effect \(d_{0i}\) on the intercept. \(V_bi\) is the intraspecific variance and \(V_bj\) is the interspecific variance.
The following models were fitted in a hierarchical Bayesian framework using the brms package (Bürkner 2017, Bürkner 2018 ) with 10,000 iterations, a warming period of 5,000 iterations and a thinning of 1/5, and 4 MCMC chains with different initial values. We obtained 1,000 estimates per parameter per MCMC chain.
The data was scaled before inference.
3.1 Paracou
Figure S3.1: Trace of the posterior estimates for the Paracou site
Figure S3.2: Density of the posterior estimates for the Paracou site
3.2 Uppangala
Figure S3.3: Trace of the posterior estimates for the Uppangala site
Figure S3.4: Density of the posterior estimates for the Uppangala site
3.3 BCI
Figure S3.5: Trace of the posterior estimates for the BCI site
Figure S3.6: Density of the posterior estimates for the BCI site
| Intercept (\(\beta_0\)) | Diameter (\(\beta_1\)) | Species variance (\(V_bj\)) | Individual variance (\(V_bi\)) | Residual variance (\(V\)) | |
|---|---|---|---|---|---|
| Paracou | |||||
| Estimate | 2.3e-03 | -2.3e-01 | 5.2e-01 | 5.5e-01 | 7.5e-01 |
| Estimation error | 2.3e-02 | 2.8e-03 | 1.7e-02 | 2.2e-03 | 5.7e-04 |
| Uppangala | |||||
| Estimate | 8.2e-02 | 1.9e-01 | 3.7e-01 | 6.6e-01 | 5.9e-01 |
| Estimation error | 4.7e-02 | 1.3e-02 | 4.4e-02 | 8.6e-03 | 1.9e-03 |
| BCI | |||||
| Estimate | 1.9e-01 | -2.2e-02 | 6.7e-01 | 4.1e-01 | 8.1e-01 |
| Estimation error | 4.7e-02 | 4.5e-03 | 3.5e-02 | 4.1e-03 | 2e-03 |
| Intercept (\(\beta_0\)) | Diameter (\(\beta_1\)) | Species variance (\(V_{bj}\)) | Individual variance (\(V_{bi}\)) | Residual variance (\(V\)) | |
|---|---|---|---|---|---|
| Paracou | |||||
| Estimate | 2.3e-03 | -2.3e-01 | 5.2e-01 | 5.5e-01 | 7.5e-01 |
| Estimation error | 2.3e-02 | 2.8e-03 | 1.7e-02 | 2.2e-03 | 5.7e-04 |
| 95% interval | -4.1e-02 - 4.8e-02 | -2.4e-01 - -2.3e-01 | 4.9e-01 - 5.6e-01 | 5.5e-01 - 5.6e-01 | 7.4e-01 - 7.5e-01 |
| R-hat | 1e+00 | 1e+00 | 1e+00 | 1e+00 | 1e+00 |
| Bulk ESS | 5.4e+02 | 2.5e+03 | 9.2e+02 | 2.6e+03 | 3.8e+03 |
| Tail ESS | 9e+02 | 3.2e+03 | 1.9e+03 | 3.1e+03 | 3.8e+03 |
| Uppangala | |||||
| Estimate | 8.2e-02 | 1.9e-01 | 3.7e-01 | 6.6e-01 | 5.9e-01 |
| Estimation error | 4.7e-02 | 1.3e-02 | 4.4e-02 | 8.6e-03 | 1.9e-03 |
| 95% interval | -1.1e-02 - 1.8e-01 | 1.6e-01 - 2.1e-01 | 2.8e-01 - 4.6e-01 | 6.4e-01 - 6.7e-01 | 5.9e-01 - 6e-01 |
| R-hat | 1e+00 | 1e+00 | 1e+00 | 1e+00 | 1e+00 |
| Bulk ESS | 1.8e+03 | 1.7e+03 | 1.9e+03 | 2.3e+03 | 4.2e+03 |
| Tail ESS | 3e+03 | 2.5e+03 | 2.9e+03 | 3.2e+03 | 3.9e+03 |
| BCI | |||||
| Estimate | 1.9e-01 | -2.2e-02 | 6.7e-01 | 4.1e-01 | 8.1e-01 |
| Estimation error | 4.7e-02 | 4.5e-03 | 3.5e-02 | 4.1e-03 | 2e-03 |
| 95% interval | 9.6e-02 - 2.8e-01 | -3.1e-02 - -1.3e-02 | 6e-01 - 7.4e-01 | 4e-01 - 4.2e-01 | 8.1e-01 - 8.2e-01 |
| R-hat | 1e+00 | 1e+00 | 1e+00 | 1e+00 | 1e+00 |
| Bulk ESS | 4.2e+02 | 3.2e+03 | 1.3e+03 | 3.2e+03 | 3.7e+03 |
| Tail ESS | 8.2e+02 | 3.2e+03 | 2e+03 | 3.7e+03 | 3.9e+03 |
Figure S3.7: Variance partitioning of individual growth in three tropical forest sites (see Table above).
Overall, a large part of variability in tree growth was imputable to individual effects in the three sites, even after accounting for the effect of diameter on tree growth, showing a high intraspecific variability in growth in these tropical forests (Figure 3.7).
4 Moran’s I analysis
In order to test if individual tree growth was spatially autocorrelated within each site, we performed a Moran’s I test on individual mean growth for each species with more than five individuals. For very abundant species (with more than 3,000 individuals), we sampled 3,000 individuals with a uniform probability. Moran’s I test function was re-written from the ape package’s Moran.I function with a one-tailed test (Paradis and Schliep 2019) in order to be able to select the pairs of neighbours that are less than 100-m apart and thus assess local spatial autocorrelation.
We then computed the proportion of species with a significantly and positively spatially autocorrelated growth (with a 0.05 alpha risk) and the proportion of species with non significant spatial autocorrelation, as well as the corresponding proportion of individuals. The proportions are computed relative to the species on which the test was performed.
| Significant | Not significant | |
|---|---|---|
| Paracou | ||
| Proportion of species | 31 | 69 |
| Proportion of individuals | 78.9 | 21.1 |
| Uppangala | ||
| Proportion of species | 18.5 | 81.5 |
| Proportion of individuals | 45.3 | 54.7 |
| BCI | ||
| Proportion of species | 20.1 | 79.9 |
| Proportion of individuals | 54.7 | 45.3 |
31%, 18%, and 20% of the species showed a positive and significant spatial autocorrelation in tree growth, in Paracou, Uppangala, and BCI respectively. (Table 4.1).
Although the three forest sites are located in very different areas across the tropics, with contrasting biogeographical history, climate, regime of perturbations and topography, they all showed a substantial spatial autocorrelation of tree growth intraspecific variability, which supports the generality of our results.The magnitude of spatial autocorrelation varied across the three sites however, with Paracou presenting the biggest proportion and Uppangala the lowest, which can be related to each site specificities. The plots of the Paracou site were subjected to disturbance treatments that created numerous large gaps (Molino and Sabatier 2001). This resulted in a high heterogeneity and spatial structuration of the environment. When restricting our analysis to undisturbed control plots (see 6.1 below), the proportion of detected spatial autocorrelation in conspecific growth drops from 79% to 55% individuals, supporting the hypothesis that the more important spatial autocorrelation in Paracou is due to the stronger spatial structuration of the micro-environment associated with disturbance treatments. Conversely, the disturbance regime in BCI is known to produce frequent small gaps (Hubbell 1999), while the mature forest in Uppangala has been reported as barely disturbed (Pélissier et al. 2011). We acknowledge that the hilly topography of Uppangala could make our estimates of distances among individuals based on pairs of tree coordinates provided by field inventory less accurate than in the other two flatter sites.
5 Is variance within conspecifics inferior to variance within heterospecifics ?
In order to test if the performance of conspecifics was locally more similar than that of heterospecifics in the three tropical forests, for each species with more than five individuals and more than five heterospecific neighbours within 100-m radii around these individuals , we computed the semivariances of growth of individuals that are less than 100-m apart, considering conspecifics only on the one hand, and heterospecifics on the other hand. In the first case, semivariance was estimated as the mean of the squared difference in individual mean growth for all pairs of conspecific individuals. In the second case, semivariance was estimated as the mean of the squared difference in individual mean growth for all pairs of individuals with an individual of the focal species and one of another species. For each species, we compared the semivariances obtained with conspecifics only to the ones obtained with heterospecifics using a Mann-Whitney test.
For each site, the proportion of species with a significantly higher or lower semivariance than the one estimated with heterospecifics or for which the difference was not significant, and the corresponding proportion of individuals were then computed.
| Intraspecific variability < interspecific variability (i) | Intraspecific variability ~ interspecific variability (ii) | Intraspecific variability > interspecific variability (iii) | |
|---|---|---|---|
| Paracou | |||
| Proportion of species | 60.7 | 40.7 | 0.667 |
| Proportion of individuals | 88.8 | 10.9 | 0.28 |
| Uppangala | |||
| Proportion of species | 42.2 | 62.2 | 4.44 |
| Proportion of individuals | 57.7 | 23.6 | 18.8 |
| BCI | |||
| Proportion of species | 49.1 | 47.8 | 3.14 |
| Proportion of individuals | 76 | 19.3 | 4.69 |
| Intraspecific variability < interspecific variability (i) | Intraspecific variability ~ interspecific variability (ii) | Intraspecific variability > interspecific variability (iii) | |
|---|---|---|---|
| Paracou | |||
| Proportion of species | 60.7 | 38.3 | 1 |
| Proportion of individuals | 86.7 | 12 | 1.26 |
| Uppangala | |||
| Proportion of species | 42.2 | 55.6 | 2.22 |
| Proportion of individuals | 64.2 | 17.4 | 18.4 |
| BCI | |||
| Proportion of species | 59.7 | 39 | 1.26 |
| Proportion of individuals | 88.7 | 10.7 | 0.516 |
Figure S5.1: Histograms of the abundances of species which present a significant (dark blue) or insignificant (light blue) Moran’s I test (upper panels) and semivariance comparison test (bottom panels). The dotted lines represent the mean abundance across species with significant (plain) and insignificant (dashed) test in each case. Using a simple generalized linear model with a logit link function, we also showed that the abundance of species is a significant predictor of test significance in both cases (McFadden’s R-squared value = 0.29 for Moran’s I and 0.17 for semivariance comparison for Paracou, 0.11 and 0.25 for Uppangala, 0.19 and 0.13 for BCI; p-value of the effect of abundance on significance = 6e-12 and 5.8e-07 for Paracou, 0.036 and 0.015 for Uppangala, 2.8e-05and 0.00016 for BCI).
In all three tropical forest sites however, growth semivariance was locally significantly higher among heterospecifics than among conspecifics in a large proportion of the tested cases (61%, 42% and 49% of the species, and 89%, 58% and 76% of the individuals in Paracou, Uppangala and BCI respectively; Table 5.1).
Note that, to control for a potential effect of species abundance on the values of semivariances with heterospecifics, we replicated the analysis by computing the semivariances with heterospecifics by sampling a maximum of ten individuals per species. The results were qualitatively unchanged (Table 5.2).
6 Other tests and verifications
6.1 Undisturbed plots in the Paracou site
One of our main hypotheses was that since many environmental variables are spatially structured and that tree growth is largely influenced by environmental variables, tree growth should be spatially structured. Our analysis using Moran’s I test showed that tree growth was indeed spatially structured at our study sites. The Paracou dataset offers the opportunity to test this hypothesis further. Indeed, some plots were disturbed in the early eighties, creating artificial gaps, whereas others were not disturbed. As the creation of gaps results in a strong spatial structure of the light available under the canopy, we hypothesized that growth should be less structured in plots that were not disturbed.
| Significant | Not significant | |
|---|---|---|
| Proportion of species | 18.1 | 81.9 |
| Proportion of individuals | 54.7 | 45.3 |
Considering the undisturbed plots only, the proportion of individuals of species with a significant Moran’s I test was 55%, much lower than when including the disturbed plots (79%, see Table 4.1; Table 6.1). This corroborates our hypothesis that the openness of the canopy triggered important spatial structure in disturbed plots due to light gaps, and as in these gaps, trees tend to grow faster, the spatial structure of growth is stronger.
6.2 Smaller stems in the BCI site
In order to be able to compare all three datasets together, we removed all stems that had a DBH inferior to 10 cm in the BCI inventories, which include all stems with DBH \(\ge\) 1 cm. We replicated our analysis on the spatial autocorrelation of tree growth in the complete dataset and found that the spatial structure in individual growth was even more significant.
7 Code implementation
The whole analysis was conducted using the R language (R Core Team 2021) in the Rstudio environment (RStudio Team 2021). The tables were made with the kableExtra package (Zhu 2021), the figures with the package ggplot2 (Wickham 2009), and the code uses other packages of the Tidyverse (Wickham et al. 2019) (tidyr (Wickham 2021), dplyr (Wickham et al. 2021), readr (Wickham and Hester 2020), lubridate (Grolemund and Wickham 2011), glue (Hester 2020), magrittr (Bache and Wickham 2020)) and other R packages (here (Müller 2020), sp Bivand et al. (2013)). Some functions were implemented in C++ thanks to the R packages Rcpp (Eddelbuettel and François 2011, Eddelbuettel 2013, Eddelbuettel and Balamuta 2018) and RcppArmadillo (Eddelbuettel and Sanderson 2014). The pdf and html documents were produced thanks to the R packages rmarkdown Xie et al. (2020), knitr Xie (2014) and bookdown (Xie 2017).